The quasi-periodic Bose-Hubbard model and localization in one-dimensional cold atomic gases 



G. Roux,''Bt. Barthel/ I. P. McCulloch,^ C. Kollath,^ U. Schollwock,' and T. Giamarchi"* 

'institute for Theoretical Physics C, RWTH Aachen University, D-52056 Aachen, Germany 
^School of Physical Sciences, The University of Queensland, Brisbane, QLD, 4072, Australia 

^Centre de Physique Theorique, Ecole Polytechnique, 91128 Palaiseau Cedex, France 
DPMC-MaNEP, University of Geneva, 24 Quai Ernest-Ansermet, CH-1211 Geneva, Switzerland 

(Dated: July 24, 2008) 

We compute the phase diagram of the one-dimensional Bose-Hubbard model with a quasi-periodic potential 
by means of the density-matrix renormalization group technique. This model describes the physics of cold atoms 
loaded in an optical lattice in the presence of a superlattice potential whose wave length is incommensurate 
with the main lattice wave length. After discussing the conditions under which the model can be realized 
experimentally, the study of the density vs. the chemical potential curves for a non- trapped system unveils the 
existence of gapped phases at incommensurate densities interpreted as incommensurate charge-density wave 
phases. Furthermore, a localization transition is known to occur above a critical value of the potential depth 
V2 in the case of free and hard-core bosons. We extend these results to soft-core bosons for which the phase 
diagrams at fixed densities display new features compared with the phase diagrams known for random box 
distribution disorder. In particular, a direct transition from the superfluid phase to the Mott insulating phase is 
found at finite V2. Evidence for reentrances of the superfluid phase upon increasing interactions is presented. We 
finally comment on different ways to probe the emergent quantum phases and most importantly, the existence of 
a critical value for the localization transition. The later feature can be investigated by looking at the expansion 
of the cloud after releasing the trap. 

PACS numbers: 03.75.Lm, 61.44.Fw, 67.85.Hj, 71.23.-k 



Disordered media are known to allow for the localization of 
waves in many physical systems, both quantum and classical. 
As demonstrated by Anderson Ilili2|], increasing disorder in- 
duces a transition to an insulating state. The occurrence of this 
Anderson transition strongly depends on the dimensionality of 
the system: in one-dimension, a localized phase is expected 
as soon as disorder is present [SJ. One of the key question 
in the field of strongly correlated systems is the interplay be- 
tween interactions and disorder Using field theoretical meth- 
ods [4, 5], it was shown that, for one dimensional systems of 
bosons and fermions, interactions can lead to a localization- 
delocalization transition. For one dimensional or higher 
dimensional |6] bosons, the combination of interactions and 
disorder leads to a transition between a superfluid phase for 
weakly repulsive bosons and a localized phase (Bose glass) 
for strong repulsion. When an additional commensurate po- 
tential is present, there is a competition between the three pos- 
sible phases, namely the superfluid (SF) phase, the Mott insu- 
lating (MI) phase, which occurs for commensurate fillings and 
large interactions, and the so-called Bose-glass (BG) phase, 
which is induced by disorder. Numerical studies ^ [s | sup- 
ported the general picture and provided phase diagrams |T3l 
in one dimension where mean-field theory fails. However, ex- 
perimental set-ups in solid state physics lack a good control 
of the interactions and the disorder strength. More recently, 
cold atomic gases offered the possibility of a fine-tuning of 
the hamiltonian parameters in particularly clean set-ups. As 
a paradigm for strongly interacting gases, the SF-MI phase 
transition was demonstrated using an optical lattice A 



fine-tuning of the disorder strength is likewise conceivable. 
In this direction, several proposals were put forward: the use 
of a laser speckle iH 0, |il HI [H O IH, the use of 
heavy atoms, which provide a quasi-static potential for lighter 
atoms 1 19, 20], and finally the addition of a superlattice poten- 
tial with a wave length incommensurate with that of the lattice 
potential I^IHIlilli. 

This paper is devoted to the study of the latter situation, 
the so-called bichromatic set-up, for which experiments have 
recently been carried out ll23l 0, 0, . The one-particle 
Schrodinger equation with an incommensurate lattice has 
been widely studied f2lll2j 



and was found to 
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exhibit anomalous diffusion properties II32I1 . The main result 
of these studies, as we will recall later on, is the existence of a 
critical value of the potential above which localization occurs. 
For the many-body physics, a weak-coupling treatment of the 
potential was carried out using bosonization [33]. Quasi- 
periodic potentials were found to have an intermediate behav- 
ior between commensurate ones, and disordered ones. Exact 
numerical results on the Bose-Hubbard model with a quasi- 
periodic potential already exist ll34l [ssi [36ll but are limited to 
small systems and thus could not investigate the nature of the 
transition nor their precise location. The physics of the Bose- 
Hubbard model with a periodic superlattice has been investi- 
gated [37ti38|,i39(] and a "weakly superfluid" phase at large po- 
tential depth was found ll39ll . Very recently, Roscilde |40] car- 
ried out a more detailed study using quantum Monte Carlo cal- 
culations and a "random atomic limit" approach. This study 
gave results on the bulk system for a special choice of parame- 
ters, and an accurate description of the physics of the trapped 
cloud, focusing on static observables. A particularly impor- 
tant point in which we go beyond Ref. 40 is the description of 
the phase diagrams of the bulk system at fixed densities which 
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FIG. 1: (Color online) Phase diagrams of the bichromatic Bose-Hubbard model for densities n = 1, r (the ratio of the potential wave lengths) 
and n = 0.5. The diagrams are shown as a function of the interaction strength U and the bichromatic potential strength V2, both normalized 
by the hopping J (lines are guides to the eyes). SF stands for the superfluid phase, MI for the Mott-insulating phase, BG for the "Bose-glass" 
phase (meaning localized but with zero one-particle gap) and ICDW for incommensurate charge-density wave phase. The U — V2 line on the 
phase diagram with n = 1 indicates the J = limit for which the gap of the one-particle excitation vanishes. Black error bars are deduced 
from calculations averaging over the phase-shift (cf. Sec. II At and finite-size scaling (see Figs. [T2] and [13] for details on the n = 1 phase 
diagram). Grey error bars are roughly evaluated from calculations on systems with L = 35 and fixed = (see Figs. [TT] and [14] for details). 
In the phase diagram with density n = 1, the darker (violet) region in the BG phase is localized but could have a small gap which cannot be 
resolved numerically. 



are essential to understand the interplay between the compet- 
ing orders at stake. 

Our main motivation is to address the shape of these fixed- 
density phase diagrams for a one-dimensional system using 
the density-matrix renormalization group (DMRG) algorithm 
(see section II CI for details), which results are interpreted 
within the framework of the Luttinger liquid theory. We focus 
on the differences and similarities between the deterministic 
bichromatic lattice potential and a truly random one, usually 
consisting of a random box distribution (RBD) and for which 
the phase diagrams without a trap are known An ac- 

count of our main results are depicted in Fig.[T] which gathers 
the phase diagrams for three typical densities as a function of 
the interaction strength U and the disorder potential strength 
V2 (see Sec. II Al and Eq. (O for a precise definition of the 
hamiltonian). n is the density of bosons and r is the ratio of 
the employed lattice wave lengths, which characterizes the in- 
commensurability of the potential. A first interesting result 
is that a finite V2 > 4 is always required to stabilize the BG 
phase. We must precise that the term BG is used to call a lo- 
calized phase which is compressible (with a zero one-particle 
gap), but the detailed features of the BG phase of the bichro- 
matic potential differ from the usual RBD BG phase as will 
be discussed in what follows. Contrary to the RBD phase di- 
agram, we argue, based on numerical evidence, that there is 
no intervening BG phase between the SF and the MI phase 
at density n = 1. An incommensurate charge-density wave 
(ICDW) phase - referred to as the incommensurate band in- 
sulator (IB I) phase in Ref.'40'- emerges at finite V2 for a den- 
sity n ~ r. Lastly, we observe that the larger the density, the 
larger the extension of the SF phase is. 

The paper is organized as follows: in section |I] we first 
give the conditions under which the hamiltonian describing 



the many-body physics simplifies into a simple lattice hamil- 
tonian used for numerical calculations. We then discuss one 
of the strongest differences compared to a random box dis- 
tribution which is the emergence of plateaus in the density- 
chemical potential curve (sectionilDl. We next discuss, in sec- 
tion ini] the competition between the disorder potential and 
the interactions by computing the phase diagrams at integer 
density one and for a density for which an ICDW plateau oc- 
curs. Lastly, section HV] is dedicated to the possible relevant 
experimental probes of localization by focusing on the out-of- 
equilibrium dynamics of the system. 

I. THE BOSE-HUBBARD MODEL WITH AN 
INCOMMENSURATE SUPERLATTICE 

A. Energy scales hierarchy: Validity of the model 

This section gives qualitative arguments on the hierarchy of 
energy scales which leads to a simple lattice hamiltonian that 
captures the physics of cold bosonic atoms experiencing two 
optical lattice potentials with wave vectors ki , k2 and ampli- 
tudes Vi, V2. Similar considerations were given recently in 
Ref. im For the sake of clarity, we keep the following discus- 
sion on the regime of experimental parameters under which 
the lattice hamiltonian under study is valid. The potential en- 
ergy in a one-dimensional setup of two standing-waves is the 
Harper potential 

V{x) = Vi cos^ikix) + V2 cos^{k2X + 0) , (1) 

which is sketched in Fig.|2] A constant phase (j) is introduced 
to shift the second lattice with respect to the other, and the 
wave vectors fci and k2 can take any value. We work in the 
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limit of a large depth Vi ^ Eri for which we can restrict 
ourselves to the lowest Bloch band {Eri — (hki)'^ /2M is the 
recoil energy associated with the first lattice and M is the mass 
of the atoms) and in a situation where one intensity is much 
larger than the other, Vi ^ V2. An exact derivation of the 
lattice parameters of the hamiltonian should resort to numeri- 
cal calculations as described in Refs. |40|. Our motivation 
is to evaluate the physical effect of the perturbing potential to 
deduce the relative magnitudes of the different terms. 

To proceed, we neglect the effect of the trap on the lo- 
cal chemical potential and displacements, meaning that we 
consider the realistic situation for the bulk physics with 
A1A2/IA1 — A2I, Ai, A2 <C y/h/AIuj with uj the trap frequency. 
If V2 — 0, the effective model is the Bose-Hubbard model lUlll 
with hopping Jo and on-site interaction Uq 

Ho = -JoJ2\-^]+ib,+b]b,+i] + UoJ2n,{n,~l)/2. (2) 

3 j 

hj is the operator that creates a boson at site j corresponding to 
the minimum of the lattice potential. The local particle num- 
ber operator reads nj = b^jbj. The dependence of the param- 
eters Jo and Uq upon Vi, Eri,ki and the scattering len gth a 
can be evaluated numerically or analytically in this limit llllll . 
We now qualitatively discuss the effect of V2 (x) to the lowest 
order in e = V2/V1. 

Perturbation of the chemical potential - First, to zeroth or- 
der in e, the minima of V{x) are located at kiXj — uj + 7r/2 
with i integer and 7r/2 can be absorbed in the redefinition of 
(j). Since Vi{xj) = + O(e^), we have: 

V{xj) = V2 cos^ {rnj + (j)) = ^yll + cos (2r7rj + 2(/))] 

The important number which characterizes this bichromatic 
potential is the ratio of the wave vectors r — fc2/fci. If r is 
a rational number p/q, the hamiltonian is q-periodic. For ?- 
irrational, it has no translational invariance and the V{xj) can 
take any value in [0, V2] in a deterministic way: the result- 
ing bounded distribution is sketched in Fig. |2] The chemical 
potential thus shares features with a one-dimensional quasi- 
crystal. The order of magnitude of the coefficient of this term 
is of course V2 — eVi. This term can be larger than Jo or 
Uq, even though e <C 1, because of the factor Vi. Taking the 
parameters of Ref. |23l one finds that e ~ 0.003-0.12 while 
V2/J0 ^ 2.6-53.3. 

Perturbation of the hopping - It is difficult to treat the 
perturbation of the hopping exactly because one needs to 
know the non-trivial shape of the perturbed Wannier func- 
tions. However, we expect the hopping to be perturbed mainly 
because of the displacement of the local minima and be- 
cause tunneling depends exponentially on the distance. We 
only consider the term associated with the perturbation of the 
minima xj and Xj+i and assume a typical exponential de- 
pendence ifTTIl for the hopping Jjj+i ~ e^'^'^J+^^^j' with 

h = -^^Z^, valid for Vi ^ Eri- The modulation induces a 
slight fluctuation Sxj at site j which, to the lowest order in e. 
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FIG. 2: (Color online) The bichromatic potential (full line) with the 
same parameters as in the experiment of Ref. 23. Dashed lines show 
the two beating standing waves from which the bichromatic potential 
originates. We observe that not only the depths of the potential wells 
fluctuate but so do their positions and their width. The energetic 
landscape displays wells of typical width 1/(1 — r). The plot on 
the right-hand side shows the bounded distribution which behaves as 
1 / y''x{V2 — x) and is thus peaked around and V2 ■ 

reads: 

6xj — ~"2^^^ (2r7rj + 20) , 

so that the distance between two neighboring sites is altered 
as: 

TT T 

Xj+i ~ ^ J, sin(7rr) cos (r7r(2j + 1) + 20) . 

Hence, to the lowest order in e, we may write for this term: 

JjJ+i = ^0 [1 + eAcos (r7r(2j + 1) + 20)] 

with, up to approximations, A = ^ sin(7rr). We 
write J2 = eAJo- Even though A could be large be- 
cause of ^Vi/Eri, the factor e ensures that J2 can be made 
much smaller than Jq. More precisely, taking the parame- 
ters of Ref. 23| and using the above approximation, one finds 
J2/J0 — 0.002-0.1, the latter occuiTing for very large V2, 
much larger than the ones used in this paper Numerical calcu- 
lations of the hamiltonian parameters i22Ll40(l confirm that the 
magnitude of J2 is small within our approximations. Hence, 
we will shorten the notation to J = Jo from now on. An- 
other feature which results from this approximation is that 
the Jjj+i have the same typical fluctuations cos(27rrj) as the 
V{xj ) which is observed numerically in Ref.l40l yet for rather 
large e. 

Perturbation of the local interaction - In the deep well 
limit, the bare interaction Uq is obtained by the relation ifTlll : 

Uo-]l-kiaEri [-^J 

where a is the scattering length. Note that the ratio Uq/Jq ^ 
exp{2y/Vi/Eri) increases exponentially with the lattice 
depth li 1 IJ . This result can be obtained by approximating the 



4 



bottom of the cosines with a parabola and using gaussian Wan- 
nier functions as the simplest approximation. The fact that Uq 
increases with Vi simply corresponds to the squeezing of the 
parabola. This squeezing may also be changed at first order 
by V2. To give a rough estimate, we can compute the second 
derivative of Eq. ([T]) and obtain for the on-site interaction: 

3 

Uo + U2 coa{2rnj + 20) with U2 = ^£^^^^0 • 

Here again, since e can be tuned to be very small, one can 
work within the U2 ^ Uq regime. The perturbation of the 
on-site interaction can thus be neglected and we will use the 
shorter notation U = Uq in the following. We also note that 
the fluctuations of the local interactions have roughly the same 
cosine dependence as the chemical potential. 

To conclude, in the deep well limit Vi ^ Eri, the follow- 
ing hierarchy of energy scales 

J2,C/2 < Uo,Jo,V2 < Vi 

can be realized experimentally. Thus, we assume the corre- 
sponding lattice model for the bichromatic optical lattice: 

^ = -^Et^Vi^.- +h-c.] + Uj2n,{n, - l)/2 
3 3 

3 

+ yEO'"Jo)'"^ (3) 

3 

with Jo — {L + l)/2 the center of the trap. In what fol- 
lows, results for the phase diagrams are computed with u) ~ 
0. The trap confinement is added in a few illustrating fig- 
ures and more importantly, for the preparation of the out-of- 
equilibrium state in the study of dynamics (Sec. lIVb . 

We can briefly comment on the distribution of the on-site 
potential energies as it is the first difference with the RBD. We 
shall use the short-hand notation for the bichromatic potential 
V2{xj) = Vj = V2[l + cos(2r7rj + 2(/))]/2. The distribution of 
the Vj with an irrational r behaves as 1 / x{V2 — x) which 
diverges close to and V2 and is symmetrical with respect to 
V2/2 (see Fig.|2]l but is relatively flat at the center. Thus, this 
distribution qualitatively lies in between a RBD and a binary 
one. Its auto-correlation function reads: 

VjV~^d -V^^ = {V2f cos{2nrd + 20)/8 

where the over-bar means averaging over all sites j. The po- 
tential is thus deterministic and correlated. Though trivial, 
this remark stresses the fact that the very features of the lo- 
calization mechanism under study originates from the quasi- 
periodicity rather than the distribution itself. For instance, an 
uncorrected disordered potential with the same distribution 
would induce localization as soon V2 is finite, which is not true 
for the bichromatic one. As sketched in Fig.|2]by black and 
dashed grey lines, wells develop over a characteristic length 
scale 1/(1 — r) ~ 4.4 sites which comes from the beating of 
the two periods 1 and r of the two lasers. 



Working with finite systems raises the question of tak- 
ing the thermodynamic al limit. The system length L is 
given in units of the first lattice spacing Ai/2. First, 
to what extent can an irrational number r be approx- 
imated by a rational number? This can be answered 
by looking at its continuous fraction decomposition ll4Tll 
which gives the successive best rational approximations. 
From Ref. 23, r = 830.7/1076.8 = 0.77145245... 
is a realistic "irrational" parameter as 8307 and 10768 
are coprimes. The successive best rational approxima- 
tions ai-e 3/4, 7/9, 10/13, 17/22, 27/35, 908/1177, . . . which 
gives the lengths L — 13, 22, 35, 1177, best "fitting" the po- 
tential for non-trapped systems. As 27/35 is akeady a fairly 
good approximation of the "irrational" r of the experiments, 
multiples of 35 such as 70, 105, can be used as well. We will 
also use other lengths L and we have checked that the physics 
does not change qualitatively if the system size does not per- 
fectly fit the potential. The fact that 35 is a rather large period 
ensures that r is not too close to a simple fraction which would 
induce strong commensurability effects on finite systems. In 
what follows, we choose to work with the experimental pa- 
rameter r — 0.77145245 as Roscilde did |40] to be as close as 
possible to the experiments but we expect the general picture 
to remain true for any irrational number Furthermore a phase 
shift (j) enters in the hamiltonian and, though we expect the po- 
tential to be self-averaging for fixed (p, averaging over </> can 
help recover the thermodynamical limit. This averaging will 
be denoted by () ^ in the following. As experimental set-ups 
generally consist in an assembly of one-dimensional cigar- 
shaped clouds with different lengths (see Fig. 3 of Ref. [U 
for instance), clouds with different lengths would effectively 
experience a different phase shift (p. Furthermore, from one 
shot to another, the tubes experience slightly different poten- 
tials. This is due to the difficulty to lock the position of the 
cloud in the trapping and optical lattice potentials over several 
shots. Consequently, cf) may fluctuate from one preparation to 
another 

The last crucial parameter in the physics of the system is 
the density which plays an important role as we will see. For a 
non-trapped system, we use the notation n = N/L with N the 
total number of bosons which is kept fixed as we work in the 
canonical ensemble. For a trapped system, the local density 
varies as one moves away from the middle of the trap and the 
thermodynamical limit is recovered for uj ^ keeping N-Ju 
constant. Roscilde i40tl gave a detailed analysis of the static 
properties in the presence of a trap. Our focus is more on the 
phase diagram of the model which is always understood to be 
in the thermodynamical limit. More details with respect to 
experimental probing will be given in section |IV] All results 
of the paper are for zero temperature. 



B. Low-energy approach: bosonization 

We briefly review known results from the low-energy ap- 
proach (close to a hydrodynamic description) which will be 
useful for the interpretation of the numerics and offer a com- 
plementary point of view on the physics. The low-energy 
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physics of interacting bosons in ID are described using Hal- 
dane's harmonic fluid approach ll42ll43ll44ll in which the den- 
sity operator is expanded as: 



1 



E 



ip(27Tnx — 2(p(x)) 



(4) 



p=0,±l,±2,... 



where n is the boson density which encompasses the lat- 
tice spacing d. The effective hamiltonian of the system 
has generically a quadratic part which includes a kinetic en- 
ergy term ^ iP, with 11 = :jV6' the conjugate field of cj) 
(the commutation relations n(x')] = id{x — x') and 

[8{x), iV(/)(a;')] — iS{x — x') hold), and a density-density 
like interaction term ^ (Vt/))^. Two Luttinger parameters u 
and K give a simple parameterization of the quadratic part of 
the hamiltonian: 

n = j ^^{uK{.Iif + ^{V^f] (5) 

where u has the dimension of a velocity and K is dimension- 
less. For free bosons, only the first term remains which would 
formally correspond to the K oo limit and u = ^ mr /M 
is the sound velocity. Taking into account a local interaction 
■jp(x)^ like in the Bose-Hubbard model, the Luttinger pa- 
rameters read u — yJrdJjM and K — i:\J n/MU in the limit 
U J. When interactions are large, higher harmonics in 
the density operator have to be taken into account to describe 
correctly the local fluctuations and not only the long-distance 
ones. In the U — oo limit, i.e. for hard-core bosons (HCB), 
one obtains K = 1 as one would find for free fermions. The 
strong interaction, i.e. the second term in Eq. Q, acts as a 
Pauli exclusion term. We thus generically have 1 < A' < cx3 
for on-site repulsive interactions. The effective hamiltonian 
(|5]l provides the general low-energy description of the SF 
phase which can undergo various instabilities. 

By bosonizing the standard Bose-Hubbard model (|2|i, com- 
mensurability effects can arise from the higher harmon- 
ics Hi: 



p{x)^ =n^ + ^(V(/))^ + E cospTrpna; - 2p<j){x)] + ... 

p>i 

^From studying the renormalization group (RG) flow equa- 
tions, it is known that cosine terms such as 



U / dx cos[2'Kpnx — 2p(j3{x)] 



(6) 



densities n = 1,2,3,... allow for a Mott phase for K below 
Kc — 2. For p — 2, charge-density wave phases can appear 
for half-i nteger densities but nearest-neighbor repulsion are 
required f44l l47ll48ll to get K < I. It is important to note that 
such cosine effective potential terms effectively arise from the 
interactions. The transition towards the charge density-wave 
(CDW) state with one atom every two sites [48] is associated 
with a spontaneous breaking of the translational symmetry. 
The other possibility to generate Mott phases is to artificially 
introduce a cosine chemical potential which directly couples 
to the density. Similarly, a CDW phase induced by an external 
potential is associated with an explicit breaking of the transla- 
tional symmetry. This latter solution is possible in cold atoms 
by adding a superlattice. 

Effect of a superlattice potential - We first consider a cosine 
potential V2{x) — V2 cos(Qa::) which has only one Fourier 
component at wave vector Q = 2Trr with r rational. The 
additional term reads 

dxV2{x)p{x) — — / dxcos[{2Trn±Q)x — 2(l){x)] + ... 

(7) 

As seen previously, such terms may induce a C-IC transition 
when increasing V2 if the condition n ± r e Z is satisfied. In 
particular, the superlattice potential can become relevant for 
the densities n = r,l — r,l + r,2 — r, .... Higher harmonics 
can be generated, as discussed in Sec. Ill Al If the potential 
term is not relevant, the Luttinger parameter K is, however, 
renormalized to a lower value by the potential as it happens 
with interactions. Such commensurate potentials have for in- 
stance been studied in the context of Mott transitions |47] and 
of magnetization plateaus |49]. The physics of cold atoms 
with induced commensurate CDW phases was studied in de- 
tail in Ref. [3^ Vidal et al. [i33i] generalized this result to irra- 
tional r. For the case of a quasi-periodic potential, the critical 
value Kc remains equal to 2 if the density approximately sat- 
isfies the relation n ± r e Z. If the density does not fulfil 
this condition but remains close enough, an insulating phase 
can be reached but for smaller critical value (for a spin-less 
fermion model, Kc — 1 was found from RG). 

Disorder with a random box distribution - From Refs. 10, 
m 0], the main result is that the potential is relevant below 
the critical value Kc = 3/2, whatever the density. The re- 
sulting Bose-glass phase has no one-particle gap but an ex- 
ponentially decaying one-particle Green's function due to lo- 
calization. The correlation length scales according to ^ ~ 
exp(— A/(V2 — V2)) where V2 is the critical value for the 
transition. 



can lock the density field and induce a commensurate- 
incommensurate transition Ii45.1 (C-IC) depending on the den- 
sity and of K. Working at fixed density and varying interac- 
tions, such a term is relevant only if the density satisfies the 
commensurability condition pn = 1,2,... and if K < Kc 
with Kc — 2/p^. The opening of the gap follows a Kosterlitz- 
Thouless 1 46] (KT) law ~ exp{-A/^/U - Uc) with A 
a constant and Uc the critical value. Working at fixed inter- 
action and varying the density, the commensurate phase is 
obtained for Kc — 1/p^. For instance, for p = 1, integer 



C. Numerical methods 

The hard-core bosons physics can be solved exactly using 
a Jordan- Wigner transformation which maps the model onto 
free fermions with boundary conditions that depend on the 
number of bosons. As the method has been widely described 
in the literature, we refer the reader to Refs. ll40i Isoll . This 
method is also used to investigate the out-of-equilibrium prop- 
erties [,50.1 of the cloud in section HV] 
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We use the DMRG algorithm lU H HI to treat the 
soft-core Bose-Hubbard model (O. For disordered systems, 
sweeping has proven to be particularly important to get reli- 
able results [10, . DMRG has also been used to investigate 
the physics of quasi-periodic electronic systems 1.55, .561. Our 
implementation is based on a matrix-product state variational 
formulation which enables us to start sweeping from any 
state. In practice, we have started from either a random or a 
classical state (where the particles are located according to the 
J = limit of the hamiltonian) contrary to the usual warm- 
up method. The algorithm works in the canonical ensemble 
(fixed number of particles N) and at zero temperature. We 
typically use from 200 to 400 kept states. The number of 
bosons allowed on-site is usually fixed to N^"'^ — 4 but re- 
sults for densities larger than one have also been checked with 
up to A^''°'* = 6. For [/ > 1 and V2 < 20, the classical dis- 
tribution of particles does not have more than 4 bosons per 
sites. 

A drawback of this variational method is the occasional ten- 
dency to get trapped in an excited (metastable) state with a 
slightly higher energy that is difficult to distinguish numeri- 
cally from the ground-state. Indeed, the usual measures of the 
convergence of DMRG, the discarded weight and the variance 
{{H — E)'^) are very small for these states. Systematic tests 
have been carried out in the U ^ 00 limit. Starting from 
the classical state improves convergence for small densities or 
close to one at large V2 as one would expect intuitively. Be- 
low V2 — 4, convergence is always good which can be related 
to the physics of the systems as the potential does not induce 
localization in this regime. In the case of soft-core bosons, 
we expect an enhancement of quantum fluctuations at finite 
U to help the particles redistribute more easily. Such equi- 
libration is rendered very difficult for HCB as for strong V2, 
local densities can be very close to one. Most of the data 
have been obtained for U < ¥2- Furthermore, relying on the 
variational principle, we can use the smallest of the two ener- 
gies obtained from starting either from the classical or a ran- 
dom state. Lastly, the coherence of the results obtained from 
observables computed independently, such as the correlation 
length and the one-particle gap (see section HIHi. supports the 
good convergence of the algorithm. 



II. DENSITY PLATEAUS: MOTT AND 
INCOMMENSURATE CHARGE DENSITY WAVE PHASES 

This section describes the relation between the density n 
and the chemical potential fi for a non-trapped cloud. The 
motivation is to find first the location of the compressible and 
incompressible phases. The chemical potential is computed 
via 

fi{N) = EoiN) - EoiN ~ 1) 

where Eo{N) is the ground-state energy with N bosons. If a 
plateau emerges in the curve, its width is directly related 
to the one-particle gap defined by 

A, = £;o(^ + 1) + Eo{N - 1) - 2Eo(N) 
= fiiN + 1) - ^liN) 



n 




-l+2r 
~_ 2-2r 

A-r 



FIG. 3: (Color online) Comparison of n{y?) for a random and a 
bichromatic potential (irrational r) for hard-core bosons. Plateaus 
open for the bichromatic potential, the main ones being at n = r and 
1 - r. 



Lastly, the compressibility of the system n = dn/d^ is eval- 
uated through its discretized expression as 

= L[Eo{N + 1) + En{N - 1) - 2Eo{N)] . (8) 

For a Luttinger liquid, the compressibility is simply related to 
the Luttinger parameter and the sound velocity: 



K 



(9) 



The compressibility naturally vanishes in a plateau phase. 



A. Plateaus in the hard-core boson limit 

Following Ref . setting J — gives insight into the J ^ 
U physics. This gives the width U — max(V,) + min(Vj) of 
the various Mott plateaus centered at /i/C/ = 0.5, 1.5, 2.5, . . .. 
This is due to the fact that in the limit J = one can re- 
order the energies by increasing values and therefore the n{p) 
curve which is the integrated density of states is simply linear 
between Mott plateaus for the random box distribution and 
I] > ma.x{Vj) — min(V, ). For a bichromatic lattice, we have 
/i = V2 sin(7rn/2). What happens when J is small but fi- 
nite? The density of states evolves smoothly with J for the 
random box distribution (see Fig. |3]l. For V2 — 0, the band- 
width which develops between the Mott plateaus has a width 
4 J and a cosine relation can be observed [7J because Mott 
sub-bands with cosine dispersion are well separated. On the 
contrary, for the bichromatic lattice intermediate plateaus ap- 
pear as soon as J is non-zero. This behavior is reminiscent of 
the situation for rational r and was discussed extensively for 
free fermions, which in our case would be equivalent to the 
HCB Hmit. 

The energy spectrum and the wave function properties have 
been widely studied in the literature [|27i [28l [Sj |30i [Slj Issll . 
It was shown that gaps open in the energy spectrum. If r is 
rational, there is a finite number of gaps. If r is irrational, 
there is an infinite number of gaps at large V2, the width of 
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FIG. 4: (Color online) Opening of plateaus in the n{fj,) curves with 
V2 for hard-core bosons (V2 is increased by steps of J/4 and is in 
units of J, as jj.). Above the critical value V2 = 4, the curves acquire 
a devil's staircase-like behavior: plateaus become dense. 



which strongly depends on p and q if one writes n — p/q and 
gets larger as V2 increases. We here recall the method usually 
followed: these gaps are studied by rn successive approxima- 
tions rm — Pm/qm of the irrational number r. For a given 
m, the potential is q^ periodic and we can use Bloch's theo- 
rem on super-cells of length g,„. The one-particle Schrodinger 
equation of the hamiltonian ^ reads: 



0. 



(10) 



Using Vj+q^ = Vj and Bloch's theorem 7/;^+^^^ = e 

the spectrum is obtained by solving the determinant of size 



Vi-E 

-J 



-J 
V2-E 

-J 



-J 



-J 



E 



= 



For qm = 2, this is the simple band folding mechanism which 
opens a gap at n = 1/2, with a doubling of the unit cell. More 
generally, at most q^ — 1 gaps appear in the spectrum made of 
(7„j bands (/c) with /c e [— 7r/(7,„, Tr/gm]. Examples of 

effective dispersion relations for the bichromatic potential can 
be found in Sec.|IV] Fig.|4]displays the opening of the plateaus 
for HCB with V2. A simple real-space interpretation can be 
given for the main plateau at n = 1 — r: it amounts to fill 
each well of size 1/(1 — r) with one particle (see Fig.|2]i. The 
plateau at n = r is simply obtained with the same argument 
with holes instead of particles. Putting two particles (holes) 
in each well can lead to plateau at densities 2(1 — r) and 1 — 
2(1 — r) = —1 + 2r. The fact that the main plateaus at n — 
r, 1 — r develop as soon as V2 is turned on is expected from 
the bosonization arguments of Sec. lIBI since K{V2 = 0) = 
1 < Kc = 2 for HCB. From a momentum space point of 
view, these opening are associated with umklapp processes 
with a momentum transfer Q which, modulo 2tt, gives back 
the conditions n = r, 1 — r. In perturbation theory, processes 
with larger momentum transfers can be obtained from Eq. (|7]i 
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FIG. 5: (Color online) n{fi) for soft-core bosons with a large interac- 
tion parameter U = 16 at small and large V2 computed on a system 
with L = 35 and a fixed phase-shift <{> — 0. Inset: curves at low 
density show the comparison between hard-core bosons (U — 00) 
and soft-core ones, which is good up to finite size effects. 



with higher order terms in ¥2- For instance, to second order, 
terms with transfers 2Q (corresponding to rt = 2 — 2r and 
— 1 + 2r) will appear if K < 1/2. Consequently, a finite V2 
is required to stabilize these plateaus (see also Fig.|9|l. As V2 
is increased, such processes break the spectrum up and make 
it point-like for the critical value ¥2—4 which is beyond this 
weak-coupling bosonization interpretation. Lastly, we note 
that this is particular to the Harper model. For the Fibonacci 
chain |33], the Fourier transform of the potential is already 
dense at small V2 ■ In our situation, the Fourier spectrum gets 
denser as V2 is increased. 



B. Plateaus for soft-core bosons 

We now consider the case of a finite interaction U. First 
of all, the hard-core boson limit is likely to give the coiTect 
qualitative behavior for large U. Indeed, at low densities, an 
interaction U slightly larger than V2 might be sufficient to re- 
cover the HCB physics as multiple occupancies are already 
strongly suppressed. Densities larger than one are allowed for 
soft-core bosons. For large U, we expect to find plateaus in 
between each Mott plateaus. One may recover the hard-core 
bosons band folding mechanism inside each Mott sub-bands 
(or at least for the lowest ones). These simple observations are 
coherent with the large U numerical data displayed in Fig.|5] 
A comparison with HCB results is provided in the inset of 
Fig. |5] which proves that [/ = 16 is sufficiently large to repro- 
duce the HCB physics within the first three Mott sub-bands. 

Fig. |6] gathers the results when U < V2, unveiling a more 
surprising behavior. As discussed previously, we expect the 
HCB behavior to account for the low-density part of the curve, 
which is actually observed through the rather large width of 
the n = 1 — r plateau. Indeed, because this plateau coiTe- 
sponds to one particle in each well, the effect of interactions 
is restricted to virtual processes. For higher densities, a large 
compressible phase is obtained, manifested by the smooth in- 
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FIG. 6: (Color online) Plateaus with a large V2 ~ U. Some of the 
largest plateaus found do not correspond to the hard-core boson limit. 
The two sizes L = 35 and 70 give an idea of the (weak) finite size 
effects. 



III. LOCALIZATION INDUCED BY INTERACTIONS OR 
DISORDER: PHASE DIAGRAMS 

We have seen that contrary to the standard random box 
situation, there is not only one phase (either the BG or the 
SF) between the Ml phases but a succession of phase transi- 
tions as the chemical potential is increased. This renders the 
usual [6J interpretation of the phase diagram in the (/i, J/t/) 
plane for a fixed ratio V2/U rather strenuous ||4^ as the suc- 
cession of phase transitions breaks it up into many compress- 
ible and incompressible pieces. Thus, we prefer to work at 
fixed density and varying the two competing parameters V2/J 
and U/ J. These phase diagrams were first sketched numer- 
ically in Refs. ll34l [ssi [36] but on very small systems and 
without a discussion of the boundaries and the nature of the 
transitions. We here provide a more precise determination, in 
particular, by using scaling over different sizes and averaging 
over (f) when necessary. We now describe more precisely the 
various observables used to sort out the phases. 



A. Observables 



crease of the density. From a phenomenological point of view, 
adding atoms fills the well minima. Since U is not too large, 
the effective potential coming from the combination of the 
interaction and the superlattice potential gets smoother and 
smoother. Consequently, the associated gain in kinetic en- 
ergy favors a compressible and actually a superfluid state as 
we will see in Sec. 1111 Hi where delocalization by increasing 
the density is discussed. In between those two regimes, the 
behavior is non-trivial. Strikingly, some plateaus existing in 
the HCB limit totally disappear, such as the n — r plateau, 
while others acquire a larger width. Having gaps whose size 
increases when interactions are reduced is something rather 
counter-intuitive. These plateaus result from the interplay of 
the potential and the interactions. A real-space picture was 
given by Roscilde [40] following a random atomic limit: con- 
sidering wells of typical size 1/(1 — r) separately, the fine 
structure of the energy levels for each number of atoms inside 
the wells depends strongly on U. V2 and also on J. Connect- 
ing wells with J allows for the computation of the integrated 
density of states which is 7i(/i). Though physically enlight- 
ening, this approach is quantitatively correct for rather small 
densities. Since the observed plateaus stem from the interplay 
of the interactions and the potential, we call the corresponding 
plateau phase an incommensurate charge density-wave phase. 
They appear to be the extension of both the Mott and the in- 
commensurate HCB phases at smaller U. Bosonization ex- 
plains, at least qualitatively, the mechanism of plateaus open- 
ing in the HCB limit by considering high order perturbative 
terms coming from Eq. which gives for instance the first 
two densities r, 1 — r and 2(1 — r), — 1 + 2r. At finite U and 
when U ^ V2, the situation is more involved as both terms 
should be treated non-perturbatively and on an equal footing. 
Predicting the observed densities at which these ICDW phases 
occur is thus beyond the perturbative approach. 



In addition to the compressibility, we need further observ- 
ables to sort out the different phases realized in the bichro- 
matic set-up. The first natural one is the superfluid density ps- 
It can be computed using twisted boundary conditions: 



-lapbc 



27rL- 



E 



ipbc 



(11) 



where the ground-state energies are computed for periodic 
(pbc) and anti-periodic (apbc) boundary conditions. With this 
definition, ps actually matches the superfluid stiffness. Other 
definitions fl^ contain the density of particles n as a prefac- 
tor. The superfluid density is zero in the BG, ICDW and Ml 
phases and finite only for the SF phases. In a Luttinger liq- 
uid, the superfluid density is directly related to the Luttinger 
parameters through 



uK 



(12) 



Combined with Eq. (|9]l, K can then be computed using K = 
^npsK. This numerical evaluation only requires the calcu- 
lation of energies. K can be independently extracted from 
correlation functions. For instance, the one-particle density- 
matrix or bosonic Green's function reads {blbj) where () in- 
dicates the expectation value in the ground-state. Following 
Ref. '59', we extract the contribution of the phase 9{x) fluc- 
tuations by dividing it by the local inhomogeneous densities 
n.: 



G{\^-J\) = 



{bib,) 



(13) 



The motivation for this renormaUzation stems from the ob- 
servation that the density-phase expression of the boson cre- 
ation operator is bi — \J p{xi)e~'^^'^^^\ and the fact that the 
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correlator which features superfluid properties in bosoniza- 
tion is (e'^'-^'^e^^^'^^J-'). For a translationally invariant model, 
both definitions only differ by a constant factor Since there 
is no translational invariance, one must likewise average cor- 
relations over all couples of points with the same distance 
X = I* — :/| to obtain a smooth behavior for this correlation. 
A typical plot is given in Fig. |20] In the case of the BG, 
ICDW or MI phases, the Green's function decays exponen- 
tially G{x) (X e^^l^ . In the Mott phase, the correlation length 
^ goes as the inverse one-particle gap ^ ^ l/A^. An effec- 
tive correlation length can also be computed on a finite system 
using S: 



E.G{x) 



(14) 



This gives a correct estimate of the correlation length for the 
localized phases in the thermodynamical limit up to a factor 
\/2. A divergence of ^ with L signals a superfluid state in 
which the asymptotic behavior of the Green's function is al- 
gebraic with an exponent controlled only by the parameter K 



G{x) cx 



1 



(15) 



This allows for the evaluation of K by using an accurate fitting 
procedure on a finite system with open boundary conditions. 
This is briefly described in AppendixlAl 

Characterizing the Bose condensation of the cloud is often 
done by looking at the condensate fraction fc- It is usually 
computed on finite clusters as the largest eigenvalue of the 
matrix pij = {b'jbi). No average over sites nor normalization 
by the local density is performed here. The largest possible 
value fc can reach is the number of bosons N. In the limit of 
HCB, quasi-condensation results in the scaling fdN) oc \/N. 
A finite fc is a feature of either the BG, the ICDW or the 
MI phase. Experimentally, time-of-flight measurements are 
related to the Fourier transform of pij , namely 



Im 



(16) 



Coherence of the quantum gas is deduced from the appearance 
of a narrow central peak n{k = 0). 




FIG. 7: (Color online) Ground-state wave function of the lattice 
model ([3} with a smooth trap but no interaction for increasing per- 
turbing potential V2 ■ There is a crossover from a gaussian wave func- 
tion (log scale) to an exponential one. Note that the maximum of the 
wave function in the presence of strong V2 is not centered at the mid- 
dle of the trap. 
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FIG. 8: (Color online) The averaged inverse participation ratio 
{/(V'o))(#. as a function of V2 for a bichromatic lattice and a random 
box distribution. There is a sharp transition at Vi = 4 in the thermo- 
dynamical limit (L = 420 for uo = G). A similar sharp transition is 
also found for the system with a smooth trap. 



B. Localization of free bosons 

We start with the simplest situation of free bosons, the 
[7 = limit, in which all bosons lie in the ground-state sin- 
gle particle wave function jV^). In Fig. 1 of Ref. 23l the 
structure of the trapped wave function is obtained from the 
Gross-Pitaevskii equation. Similar results are found here for 
the lattice model ^ as shown in Fig. |7] which displays the 
qualitative change of shape from a gaussian to an exponential 
structure. In order to quantify the localization transition of 
a single-particle wave function It/j), one can use the inverse 
participation ratio which is usually defined as 

3 



\j) is the state at site j in the real space basis. In the thermody- 
namical limit, /(V') goes to zero for a delocalized state with 
a typical scaling 1/L or ^/uj for respectively a non-trapped 
and a trapped system, while it remains constant for a local- 
ized wave function. Based on an exact duality transforma- 
tion of the one-particle Schrodinger equation ( fTOl i. the local- 
ization of the wave function has been conjectured by Aubry 
and Andre to happen at the critical value V2 = 4. This 
conjecture is illustrated in Fig.[8]which displays {I{iIjo))cI> as 
a function of ¥2- In comparison with the RBD evolution, the 
bichromatic set-up displays a sharp transition even for a finite 
trap frequency provided it is small enough. For the RBD, lo- 
calization occurs as soon as V2 7^ |3] with a typical scaling 
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Luttinger exponent K for liard-core Ijosons 
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FIG. 9: (Color online) The Luttinger parameter K for hard-core 
bosons as a function of V2 and the density n on a. finite system with 
L = 175 for an irrational r = 0.7714. . .. Strictly speaking, K 
should be equal to 1 in all superfluid phases and in the gapped 
phases. Up to finite size effects, vertical bands reveal the successive 
openings of gaps as V2 is increased. 



C. Localization of hard-core bosons 



We have seen that plateaus emerge in the ri(/i) curve as 
soon as V2 is turned on and that the spectrum is point-like 
above V2 — 4. The extension of the wave functions is related 
to the nature of the energy spectrum and it was shown l28ll 
that all wave-functions are extended below V2 while they are 
all localized above. Consequently, we expect the HCB to lo- 
calize above V2, whatever the density. Below V2, HCB can 
be either in a SF or in an ICDW state. To illustrate this situa- 
tion, we plot in Fig. |9] the behavior of the Luttinger exponent 
_ftr as a function of V2 and the density n. It nicely shows that 
K — lin superfluid phases as expected for HCB but vanishes 
(up to finite size effects) for the densities corresponding to the 
ICDW phases, the main ones being located at n = r and 1 — r. 
Many gaps develop as the critical point is approached and the 
shrinking of the bands renders the low-energy approximation 
and calculation of K difficult close to this point. 

As a partial conclusion, the two limiting cases U = 
and ?7 = 00 of the {U, V2) phase diagrams of Fig. [T] can 
be summed up as follow: (i) for a generic density n (mean- 
ing that it does not correspond to an ICDW plateau) and also 
for the U ^ limit whatever the density, the system re- 
mains superfluid for V2 < 4 and localizes in a BG phase for 
V2 > 4 with a correlation leng th which behaves according 
to - ln(V"2/4) miMIlM, (ii) for a density close to a 
plateau phase (for instance n = r or 1 — r) and U — 00, there 
is a transition towards an ICDW phase for a critical value of 
V2 which is smaller than 4 (and precisely equal to for n = r 
or 1 — r), (iii) for the commensurate integer density n = 1 and 
U — 00, the system remains in the MI phase ground-state for 
any ¥2- 




FIG. 10: (Color online) Superfluid-Bose glass localization transi- 
tion. The parameters U — 2 and n = 0.63 are chosen such that the 
system is in the superfluid phase at V2 =0. (a) The inverse corre- 
lation length scales roughly as \JV2~ where V2 — 6.9 is much 
larger than the non-interacting critical value V2 = 4. Below V2, we 
have the scaling ^ ~ L typical of a superfluid state, (b) ps provides 
an independent determination of the transition, (c) Evolution of the 
Luttinger exponent K as a function of V2 ■ The dashed line displays 
the RBD result Kc — 3/2. (d) Averaged one-particle density-matrix 
G{x) for increasing V2 showing the transition from an algebraic to 
an exponential decay. 



D. The superfluid - Bose glass transition for soft-core bosons 

We address here the direct transition from the SF to the BG 
phase which occurs for a generic density by increasing the 
strength of the potential V2. Fig. [TO] provides the evolution 
of the Green's function G{x) showing the localization tran- 
sition. First, a finite "disorder" strength with a critical value 
V2 — 6.9 is necessary to obtain exponentially decaying corre- 
lations. This value is larger than the U = and U — 00 lim- 
its; interactions have a delocalization effect on the BG phase 
similarly to the RBD box results. Qualitatively, this can be 
understood by starting from the localized state. There, the 
condensate is fragmented into pieces. Repulsive interactions 
will make the condensate fragments inflate and, by doing so, 
will help make them overlap and build coherence. For bosons, 
interactions thus helps delocalization. Interestingly, comput- 
ing the Luttinger exponent from the coiTelations shows that 
the critical value Kc at the transition is smaller than the RBD 
result 3/2. The scaling properties of the transition thus dif- 
fer from the standard SF-BG transition. Finding a Kc smaller 
than the RBD result for the Harper potential is well compati- 
ble with the anal ytic al finding for K found for the Fibonacci 
potential in Refs. ^y. 

To proceed with the discussion of the competition between 
interactions and the disordered potential, we compute with 
DMRG the phase diagrams of the system in the three generic 
cases n = 1 (competition between the SF, Mott and BG 
phases), n ~ r (competition between the SF, ICDW and BG 
phases), and lastly n — 0.5 (competition between the SF and 
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FIG. 11: (Color online) Phase diagram for n — 1. Observables are computed with DMRG for a system with L — 35 and fixed phase-shift 
<jf> = as a function of U and V2. The V2 = line shows the Mott transition at [/c = 3.3 while the U = line shows the free boson 
localization transition around V2 = 4. The Mott insulating phase gets qualitatively delocalized as V2 increases for U not too large. Increasing 
U delocalizes the BG phase if V2 is not too large. 



BG phases only). The summary of the phase diagrams is given 
in Fig. [U 



E. Phase diagram at n = 1 

All observables relevant for the construction of the phase 
diagram as a function of the interaction U and the potential 
depth V2 are reported in Fig. [TT] for the integer density ?i = 1. 
The Mott phase is characterized by a finite gap Ac ~ 1 a 
zero SF density and a finite and small (of order unity) conden- 
sate fraction. It emerges at the bottom right corner above the 
critical value [48] Uc ~ 3.3 for V2 = 0. We observe that Uc 
increases with V2 as for the RBD, meaning that V2 destabi- 
lizes the Mott phase. One can understand from a simple local 
on-site energies argument: the disorder reduces the minimum 
one-particle energy gap in the atomic limit. The BG phase 
exhibits exponentially decaying coiTelations, a zero SF den- 
sity and a non-diverging condensate fraction but no gap. It 
emerges for large V2 region of the V2 > U half-part of the 
phase diagram. Note that the BG has a condensate fraction fc 
that is slightly larger than for the MI phase, qualitatively due 
to the fact that coherence should remain significant over the 
typical length scale of the wells, namely 1/(1 — r). The SF 
phase has a finite SF density but no gap and algebraic coiTe- 
lations. It generically emerges at low U and low V2 and sur- 
prisingly extends into a hand print like pattern. A very small 
(compared to U and V2) one-particle gap is observed for large 
U and large V2 but we cannot conclude whether it is a finite 



size effect or not (see Fig.flj]). 

Superfluid-Mott transition and intervening Bose-glass 
phase - An important question is whether the SF and Mott 
touch each other at small but finite V2. In other words, is there 
always an intervening BG phase between the SF and the MI 
phases as for the RBD (|9l llOtl ? For the bichromatic poten- 
tial, we however have reasons to think that small V2 might 
not be as relevant as for true disorder since a large critical 
value exists for both hard-core and free bosons. To address 
this issue numerically, we have compared the scaling of the 
most relevant observables for the known case V2 — and for 
V2 = 2 (see Fig.[l2]). When V2 = 0, the SF-MI transition is 
of the Kosterlitz-Thouless type leading to an opening of the 
one-particle gap Ac cx exp{—A/^/U — Uc) above the crit- 
ical value Uc, with A a constant. Such an opening gives a 
good fit to the extrapolated data (see Fig. [T2] | but does not 
precisely give Uc- Finding Uc is rather achieved by using 
the weak-coupling RG result Kc — 2 for the KT transition. 
Fig. [T2I shows that Uc — 3.3 ±0.1 for V2 = in agreement 
with results of Ref. 48. Within eiTor bars, the scalings of the 
superfluid stiffness p,, and correlation length also agree with 
this critical point. Note that because of the very slow open- 
ing of the one-particle gap in a KT transition, the correlation 
length and superfluid density show much smoother finite size 
effects than for the SF-BG transition illustrated in Fig.[TO] For 
V2 = 2, if a BG is present in between the SF and MI phase, the 
one-particle gap Ac should open after the superfluid stiffness 
scales to zero. Up to numerical accuracy, data are consistent 
with a direct SF-MI transition of the KT type with a slightly 
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FIG. 12: (Color online) Superfluid-Mott insulator transition for n = 1. Cuts along the [/-axis for V2 ~ and 2. Left: data giving similar 
results to those of Ref . i48ii (the vertical bar being Uc — 3.33 ± 0.1). Data are computed with A'^''"'' = 4 and the Luttinger exponent K is 
determined using either ps and k or G(x) (see Sec. lIII At . Scaling of ps, £,/ L and the criteria Kc = 2 gives the same critical point within 
error bars. The scaling of the condensate fraction /c (not shown) is not simple at the transition and the HCB scaling /c oc \/iV does not hold. 
Right: the same observables for V2 ~ 2 and A'^''"'' — 3 (this cutoff induces a non-physical decrease of the superfluidity at small U, but does 
not affect the transition as seen, for instance, from the behavior of the Luttinger parameter K for A'^''"'' — 4). Insets show scaling behavior for 
ps and ^/L after averaging over several different (j>. A critical point Uc — 3.6 ± 0.1 is found which corresponds to Kc — 2. Furthermore, 
the one-particle gap (not averaged over (f)) is best fitted by a Kosterlitz-Thouless opening (we fixed Uc ~ 3.6 for the fit). Up to numerical 
precision, we infer from these results that there is a direct transition between the SF and the MI phases, with no intervening BG phase. 



larger critical interaction Uc — 3.6 ± 0.1. We observed that 
averaging over (j> is needed to ensure a good crossing of the 
scaling curves (see insets of Fig. [12]). Note that for the RBD 
situation, V2 — 2 would correspond to a disorder amplitude 
A = 1 in Ref. 10 (or A = 0.5 in Ref. 9) for which the BG 
phase already has a significant width. Interestingly, the n = 1 
phase diagram has a similar sha pe a s the one with a commen- 
surate potential with r = 1 /2 13911 . In this case, there is at 
large V2 a charge density wave phase is a gapped with two 
particles each two sites for which Kc — 1/2. A direct SF-MI 
transition is found because the term (|7j will not be relevant 
for small V2. Our results suggest that in the incommensurate 
case, the potential remains irrelevant as well. 

Superfluid-Bose Glass transition - We now turn to the dis- 
cussion of the contour of the SF-BG transition which displays 
a "hand print" pattern. First, contrary to the RBD, the BG 
phase emerges only above = 4 and for much larger values 
for small U. Secondly, V2 increases with U at small U which 
is similar to the delocalization by interactions observed in the 
RBD case. Similarly to what was found in Fig. [10] the inverse 
correlation length has a power-law behavior above the critical 
point with a Luttinger parameter smaller than 3/2. The con- 
vexity of the SF phase contours changes contrary to the RBD 
phase diagram, leading to this hand print pattern. To under- 
stand if these reentrances of the SF phase inside the BG phase 
are not a finite size effect and remain after averaging over </>, 
we show the averaged (/Os)^ and {S^)^, for various system sizes 
in Fig. [13] The behavior of ^ and ps suggests two reentrances 
of the SF phase and in particular, a sharp but clear one close to 
the transition to the MI phase. The U — V2 line corresponds to 




FIG. 13: (Color online) Reentrances of the SF phase with increasing 
interactions. Cut at V2 ~ 10 in the phase diagram of Fig.[TT] Error 
bars are related to the average over the phase-shift (f} (about twenty 
samples). The scalings of the correlation length (a) and of the su- 
perfluid density (b) suggest two reentrances of the SF phase as U is 
increased. In between, a localized phase is found and a MI phase is 
obtained at large U according to the gap (c). In the intermediate BG 
phase (for U ~ 6-9), a small gap is found, but for the chosen sizes 
L — 22, 35 and 70, the scaling ^ ~ 1/A shown in (d) seems to be 
satisfied only in the MI phase, meaning that strong finite size effects 
are present. 



the "transition" between the MI and BG phases in the atomic 
limit. It gives a rough estimate for the extension of a SF phase 
at large U and V2 which does not occur for the RBD situation. 
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FIG. 14: (Color online) Phase diagram for n ~ r. Observables are computed on a system with fixed size L = 35 and fixed phase-shift (p = 0. 
A large ICDW phase emerges at large U. A finite V2 is required to stabilize it. The SF phase extends along the V2 = contrary to the phase 
diagram with n = 1. 



We expect that the SF phase vanishes for large U and that 
there is a direct MI-BG transition around the U ^ V2 line. A 
similar emergence of the superfluid phase around the atomic 
limit was found in Ref.H^for the case of a commensurate po- 
tential where the SF phase competes a CDW and a MI phase. 
In Fig. [13] the intermediate localized phase between the two 
SF reentrances displays a small gap. This phase could have a 
finite gap but we observe that if so, it cannot be distinguished 
from finite size effects. Comparing again with the commensu- 
rate case |39], the main difference (apart from the finite gap) 
is the extension of the SF phase along the U = line. For 
the Harper model, there is no such extension because of the 
localization of the single-particle wave-function. 



F. Phase diagram close to the density n — r 

A density which satisfies the criteria n r allows for the 
realization of an ICDW phase which competes with the SF 
and BG phases. The ([/, V2) map of the observables is given 
in Fig. [14] The ICDW phase has a finite gap and exponen- 
tially decaying correlations as the MI phase. Similar qualita- 
tive features are found with the ICDW phase replacing the MI 
phase. However, a finite V2 is of course required to stabilize 
the ICDW phase contrary to the MI phase. Secondly, a finite 
V2 is needed to stabilize the BG phase. As a consequence, the 
SF phase extends to large U close to the V2 — line. As dis- 
cussed in sectionini the ICDW is a new feature compared with 
the RBD phase diagram given in Refs.loifTolfor n = 0.5. Sim- 
ilar reentrances of the SF phase into the BG phase are found 
at fixed V2 and increasing U . The U — oo line of the phase 



diagram would give an ICDW phase everywhere except for 
V2 = since the n ^ r plateau occurs as soon as V2 is finite 
in the HCB hmit. 



G. Phase diagram for a generic density n = 1/2 



Lastly, the phase diagram for a generic density n = 1/2 
has been computed to discuss only the competition between 
the SF and BG phases (data not shown, see phase diagram 
in Fig. [T]for results). We must note that, ICDW plateaus can 
however appear for generic density in a region where [/ ^ V2, 
as we show in Fig. |6]for the particular choice of parameters 
V2 — 8 and U — 6. In this case, the ICDW phase would 
have a finite domain in the (f7, V2) map (contrary to the n — r 
phase diagram), because the ICDW phase is not reaUzed in 
the HCB limit. The observables suggest that Ac remains zero 
in the whole parameter range, while the BG is bounded by the 
V2 = 4 line and the SF phase slightly extends inside the BG 
phase for small U . However, critical values V2 for the SF-BG 
transition are found to be smaller than for n ~ r, themselves 
smaller than for n = 1. The same qualitative argument stating 
that the lower the density, the closer the physics is to the HCB 
can be put forward. The SF region extends with the density 
of the system. This observation will be now more precisely 
discussed. 
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FIG. 15: (Color online) Delocalization by increasing the number of particles. Observables for increasing density when V2 — U (same 
parameters as in Fig.[6]l. There is a delocalization transition with increasing density. A few ICDW phases can be seen at intermediate fillings 
(dashed brown vertical lines in the rh.s. figure showing the one=particle gap Ac). Grey areas denote the localized regions (either BG or 
ICDW). The weakening of the superfluidity at large density might be an artefact of the cut-off in the number of bosons kept per site which is 
here iV''"' = 4. 



H. Delocalization via increasing the density of bosons 

An complementary approach to these ([/, V2) phase dia- 
grams at fixed density is to keep U and V2 constant and to look 
at the observables as a function of the density n. From dipole 
oscillations measurements. Lye et al. ll23ll observed a delocal- 
ization transition by increasing the number of particles. We 
now address the non-trivial case ofU ~ V2 by setting [/ = 6 
and V2 = 8 corresponding to the parameters of the curve 
of Fig.|6] Results for the same observables as for the phase di- 
agrams are plotted in Fig. [15] We found transitions between 
the three different phases BG, ICDW and SF. At low densities, 
double occupation for bosons is strongly suppressed because 
of the finite U . Consequently, the behavior is qualitatively 
the one HCB would have: V2 being larger than 4, localiza- 
tion exists at low densities. The superfluid density, correlation 
length and one-particle gap confirm the presence of the BG 
phase. At large densities, a SF emerges which is something 
well-known without disorder because the lobes of the Mott 
phases shrink at large densities ll48ll . In addition, the disor- 
dered potential has a tendency to reduce the size of the Mott 
phases as we have seen. Very qualitatively, some particles 
fill the wells of the disorder potential so that the remaining 
ones feel a smoother effective potential allowing for a gain in 
kinetic energy leading to superfluidity. This behavior for an 
irrational r is qualitatively similar to what was observed for 
a rational r (see Fig. 23 of Ref. |39) except that no BG, but a 
"weakly superfluid" phase is realized in this latter case. Be- 
sides this sharp BG-SF transition, peaks in the one-particle 
gap Ac uncover the presence of ICDW phases within both the 
BG and the SF phases. These phases naturally correspond to 
the plateaus in Fig. |6l 

IV. PROBING THE BOSE-GLASS PHASE WITH 
OUT-OF-EQUILIBRIUM DYNAMICS 

A. Static observables 

The question of probing experimentally the BG phase with 
respect to the other possible phases is particularly important. 



First, the simplest observable obtained after time-of-flight 
measurements is related to the momentum distribution of the 
atoms n{k). This measure helps distinguish between coherent 
and incoherent phases by looking, in particular, at the fc = 
peak. A sharp and high peak is the signature of a coherent 
phase, the superfluid phase. Because of short-range correla- 
tions, both the MI and the BG phases will give a much smaller 
peak broadened with a typical width of Fig. [T6l displays 
n{k) in the BG at small U . In addition to the central peak, 
satellites peaks at A; = ±27r(l — r) emerge as a signature of 
the underlying superlattice. However, in the experiment pic- 
tures, the Wannier envelope and the broadening of the peaks 
due to scattering events during the time-of-flight will change 
the observed shape. It is expected that the additional satel- 
lite peaks are too small to be experimentally resolved and are 




k 

FIG. 16: Typical n{k) in the Bose-glass phase of Fig.[To]in a system 
without a trap. Satellites peaks at wave vector 27r(l — r) are visi- 
ble if the disorder is not too strong. The width of the middle peak 
typically gives the inverse correlation length (a Lorentzian of width 
^ computed from Eq. il4t is given in dashed lines for a qualitative 
comparison). 
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FIG. 17: (Color online) Expansion of HCB condensates. At t = 0, 
the system is in the ground-state of the hamiltonian with the trap plus 
"disorder" potentials using a fixed chemical potential jj, — (about 
51 particles) and a trap frequency lj = 0.03. At t > 0, the trap 
confinement is switched off abruptly. Figures show the evolution 
of the local density profile as a function of time. An expansion is 
observed for systems with V2 ~ (without disorder) and V2 = 2, 
but not for a random potential and a bichromatic potential with V2 = 
6 > V'f. For the RBD, V2 defines the width of the box distribution. 
When there is an expansion, reflections on the boundaries of the box 
in which the condensate expands can be seen. 



washed out if either V2 and/or U are too large. Thus, n{k) can 
only be used to distinguish the superfluid from the Bose-glass 
or Mott insulating phases. However, it would not help distin- 
guish the MI from the BG phase. Refs. 34, 40 found a similar 
behavior and, in the second reference, a non-monotonic evolu- 
tion of the central peak n{k — 0) with increasing V2 has been 
established. The reinforcement of the superfluidity upon in- 
creasing V2 at fixed J7 in a trapped cloud must be reminiscent 
of the MI-SF-BG transitions of the phase diagrams of Fig. [T] 
Noise correlations ll60ll were proposed ll40l l6ll| as a possible 
probe for the BG phase and measured in Ref. To. However, 
this observable catches the fact that density coiTelations reveal 
the presence of the underlying superlattice ll40[|6lll but not the 
gapless nature of the excitations |40]. It is therefore necessary 
to look for additional evidence of localization. 



B. Expansion in tlie lattice potential 

As often done in experiments fl^ [ll [Tl [Tl [TtI [Tsl B, 
transport measurements are a better fashion to probe localiza- 
tion. In order to show the existence of a critical point for the 
localization, we propose to look at the expansion of the cloud 
when the trap is released |[il[ll[il[il|53,|6l|6l. Observ- 
ing the expansion in the optical lattice is a particularly appeal- 
ing experiment as the hamiltonian governing the dynamics is 
the one of the bulk system (with a; = 0) for which we have 
computed the equilibrium phase diagram. The confinement is 
used here to prepare an out-of-equilibrium state for this hamil- 
tonian. 

For the sake of clarity, we first discuss the expansion of 
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FIG. 18: (Color online) Dynamical critical point - Width of the 

atomic distribution W = \J Uj {j — jo)^ of a HCB condensate as 
a function of V2 for several increasing times (in units if the hopping). 
The width Wo at t = has been subtracted for clarity. Parameters 
are the same as in Fig. [17] Up to finite trap frequency effects (see 
text), the dynamical critical point is identical to the equilibrium one 
(at least for HCB). 



HCB. The spreading has been studied before in the HCB limit 
for homogeneous lattices ifsoll . and for soft-core bosons in 
commensurate lattices ll64ll . Fig. [17] displays the expansion 
for free HCB (V2 = 0), for two bichromatic potential ampli- 
tudes, below (V2 = 2) and above (V2 — 6) the equilibrium 
critical point, and also a situation with a RBD potential of an 
amplitude V2 = 2. The system is prepared in the ground-state 
of the hamiltonian with the confining potential (we chose a 
trap frequency of oj — 0.03). Alt — {), the trap frequency is 
set to zero and the condensate is free to expand into the lattice. 
For V2 = 0, the expansion of the edges of the condensate is 
roughly linear, with a typical velocity 2 J corresponding to the 
maximum group velocity (see below). For the RBD potential, 
the expansion is inhibited for the amplitude V2 = 2. How- 
ever, for the bichromatic set-up, the same potential strength 
does not prevent the condensate from expanding. Still, V2 — & 
induces a localization of the condensate similar to the one ob- 
served for the RBD potential. 

One may ask whether the critical point of the dynamical 
localization observed in Fig. [T7] is the same as the equilib- 
rium one. In the case of HCB, we know that all single-particle 
wave-functions are localized above V2 = 4 (see Sec. Ill Al and 
references therein). Consequently, we expect the dynamical 
critical point to be identical to the equilibrium one. To support 
this statement, we show in Fig.[T8]the width of the atomic dis- 
tribution of the condensate after several times of expansion, 
as a function of the "disorder strength" V2. The dynamical 
critical point is found to be very close to V2 = 4, within a 
±0.1 window (grey rectangle in Fig.[T8]l. One observes that, 
slightly above V2 = 4, the condensate still spreads a little bit 
with time. This may be understood as a finite trap frequency 
effect. Indeed, the edges of the condensate can spread over a 
few sites if the initial trap is too steep. In this case, the start- 
ing atomic distribution is too far from the local atomic distri- 
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FIG. 19: (Color online) Effective dispersion relations for HCB - the Fourier transform of the single-particle wave functions |^fe|^ is plotted 
as a function of the single-particle energy and the pseudo-momentum k. At the right each dispersion curve is the initial occupation numbers 
(from the state prepared using the trap confinement) v.s. the same single-particle energies for two different chemical potential /i = —1.5, 0. 
This shows which states participate to the expansion. Other parameters are the same as for the expansions shown in Fig.il7i 



bution tliat is expected (locally) in the bulk of a non-trapped 
system, and particles have to be redistributed. In the limit of 
vanishing initial trap frequency, we expect the transition to 
be sharper Note that, as the localization of all single-particle 
wave-functions in the spectrum does not depend on r (pro- 
vided it is irrational), we expect the equality of the dynam- 
ical and equilibrium critical points to hold independently of 
r. Fig. [18] may be interpreted as the vanishing of all effec- 
tive group velocities at the critical point (see below). Thus, 
the exact critical value VJf = 4 for the localization could 
be probed experimentally with this technique for free or hard 
core bosons. 

Following Ref. 21, a more precise description of the HCB 
expansion can be carried out by looking at the one-particle ef- 
fective dispersion £{k) for HCB. Without translational sym- 
metry, wave vectors k are not good quantum numbers but 
looking at the Fourier transform of the one-particle wave func- 
tions ijjk ^ J2j e^^-'i^j and plotting iV'fcP as a function of the 
pseudo-momentum k provides an effective dispersion. The 
features of the expansion depend mainly on two properties. 
First, the group velocities Vg{k) = de{k)/dk derived from 
the effective dispersion relation convey the typical maximal 
speed at which expansion evolves. Second, the expansion also 
strongly depends on the initial occupation numbers [t ~ 0) 
of the eigenstates of the hamiltonian with lo — Q. This occupa- 
tion is plotted together with the dispersion relation as a func- 
tion of the "single-particle energy" e in Fig.[T9lcorresponding 
to the expansion observed in Fig. [TT] For V2 = 2, for which 
there is no localization, the effective relation dispersion dis- 
plays gaps as we have seen from section and well-defined 
bands with a shorter periodicity originating from the band 
foldings induced by the potential (see Sec. Compared 
with the single cosine dispersion obtained without disorder. 



several shifted bands exist due to Bragg scattering with the 
potential. One can convince oneself that opening gaps low- 
ers the maximum possible group velocity. Thus, compared 
to a system with no disorder (Fig. [19] V2 = 0), the expan- 
sion for the bichromatic potential below V2 = 4 will always 
be slower if the e ~ i) state (associated with the maximum 
group velocity 2 J) is occupied in the initial state without dis- 
order This explains the qualitative features of the situations 
for which the condensate expands in Fig. [17] When V2 = 0, 
the expansion is slower when the chemical potential is below 
much (not shown), as can be guessed from the initial oc- 
cupation numbers in Fig. [19] When V2 = 2, the structure of 
the expansion is rather homogeneous at low chemical poten- 
tial (not shown) but becomes inhomogeneous and faster for 
larger ^ (see Fig. [TT}. The presence of two different speeds 
might stem from populating bands with different maximum 
group velocities as can be seen in Fig. [19] For the bichromatic 
potential with V2 = 6, no bands can be distinguished as the 
signal does not show well-defined pseudo-momentums. The 
RBD potential displays a very different effective dispersion 
faded by the disorder, but which still retains the whole feature 
of the cosine dispersion without disorder These two pictures 
illustrate that the localization mechanism for the bichromatic 
and a RBD potential is qualitatively different: the first one 
is rather associated with a band folding mechanism while the 
second rather coiTesponds to strongly scattered single-particle 
states. In this respect, one can view the "weakly superfluid" 
phasejfound for commensurate superlattices with a large V2 in 
Ref. |39| as a precursor of the Bose-glass phase of incommen- 
surate lattices. 

How these results can carry to investigate the physics at fi- 
nite U is an important and challenging question that needs fur- 
ther investigations going clearly beyond the goal of this paper 
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Indeed, from the numerical point of view, the expansion of 
strongly correlated soft-core bosons is accessible with time- 
dependent DMRG only until times of order 10 J | (64l] . while 
Fig. [TS] shows that a reasonable determination of the out-of- 
equilibrium critical point requires times of, at least, lOOJ. 
Thus, the question of the dynamical localization at finite U 
of the model ([3]) and its relation with the equilibrium phase 
diagrams of Fig. [T]remains an open question. First, the most 
naive prediction would be to expect a similar physics than the 
HCB one to hold, at least for large enough interations. Fur- 
ther qualitative arguments can be given on the expansion for 
intermediate U, for which we can use two results from the 
equilibrium phase diagrams studied in Sec. |llll (i) the criti- 
cal values to observe localization are all larger than V2 ~ 4, 
whatever U or the density, (ii) at small densities, the physics 
is essentially equivalent to the one of HCB. First consider a 
situation where V2 < 4. The starting trapped state is expected 
to have regions which can be locally SF or Ml but not EG, 
since there is no intervening BG phase. Provided V2 < 4, 
an expansion is then expected systematically (whatever U or 
the total number of particles), because the edges of the con- 
densate would be in a SF state (see examples of expansions of 
strongly-correlated soft-core bosons in Ref. |64]). ForV2 > 4, 
the situation is more subtle as, for a given U and V2, the oc- 
curence of localization depends on density in the non-trapped 
condensate (see for instance Fig.fTSl). The starting state struc- 
ture is complex and local-density approximation not necessar- 
ily valid 140(1 ■ Very qualitatively, the edges of the condensate 
would be in a localized state while, if the density at the center 
of the cloud is large enough, SF or Mott regions could also ap- 
pear. However, if expansion there is, the local density will de- 
crease with time. When the density becomes small enough to 
neglect interactions, localization would then be expected since 
V2 > 4 and one enters the HCB regime. A possible scenario 
could thus be a systematic localization after either a transient 
regime with expansion, or no transient regime. However, as- 
certaining whether the above qualitative arguments could be 
spoiled by other effects is difficult. For instance, it is known 
in the completely different limit of very weakly interacting 
bosons (Gross-Pitaevskii limit) that non-linear effects due to 
interactions can lead to some kind of localization even in a 
purely effective periodic potential i65i l66tl . How to go from 
such a limit to the relevant one for the Anderson localization 
in strongly interacting one dimensional systems is clearly a 
question that will need further experimental and theoretical 
work. 
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FIG. 20: (Color online) A typical example of fitting the averaged 
data from DMRG with conformal field theory results to extract the 
Luttinger exponent K. System size is L = 128. 



soft-core bosons, leading to a finite critical value of the lo- 
calization transition V2 > 4. The critical values found are 
high, possibly sufficiently high to allow for an experimental 
demonstration of a localization transition. In this perspective, 
static observables give a clear evidence to distinguish between 
coherent and localized phases, but their ability to sort the BG 
from the (small-F2) MI phase is less obvious. On the con- 
trary, the expansion of the condensate after switching off the 
confinement is proposed to provide a simple and rather clear 
signal to detect the localization transition. This was shown 
explicitly in the case of hard-core bosons but remains an open 
question for soft-core bosons. 
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V. CONCLUSION 

The Bose-Hubbard model with a quasi-periodic potential 
was shown to display a rich phase diagram including a Bose- 
glass phase (localized but compressible), and incommensu- 
rate charge-density wave phases in addition to the superfluid 
and Mott phases. While localization induced by this random- 
like potential is found, the underlying mechanism differs from 
the RBD situation: the band folding mechanism known pre- 
viously for free and hard-core bosons (or fermions) holds for 



APPENDIX A: METHOD TO FIT THE BOSONIC GREEN'S 
FUNCTION ON FINITE SYSTEMS 

We use conformal field theory results |43] for a system of 
length L with open boundary conditions to fit the bosonic 
Green's function defined in Eq. ( fTSl l. In the case of corre- 
lations of the type (^q''^(^)q-'^s(x )^ j^^^ . 



G{x,x') = A 



^d{2x\2L)d{2x'\2L) 
d{x + x'\2L)d{x-x'\2L) 



1 

2K 



(Al) 
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with K the Luttinger parameter, A a constant and d the con- 
formal length 



d{x\L) = - 

TT 



sm 



/Tra; 



Because there is no translational invariance, the correlations 
depend on both x and x' . Hence, in order to perform a fit, one 



has to average Eq. ( lAll i over the results obtained with fixed 
distance x' — x. Strictly speaking, formula dAlb is valid for 
1 ^ x^xWx' — x\ <C L and we have to remove the corre- 
sponding contributions. Practically, fits are rather good up to 
distances comparable with L as one can see in Fig. |20] and 
significantly improves the determination of K compared with 
a simple algebraic fit. 
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